package de.latzko.pattern.monteCarlo

import java.util.Random
import java.util.Random._
import de.latzko.pattern.math._
import de.latzko.pattern.math.NearestNeighbour
import scalala.tensor.dense.DenseVector._
import scalala.tensor.dense.DenseVector
import scalala.Scalala._
import scalala.Scalala

class Wolff(val vectors: Array[DenseVector], val q: Int, val k: Int) extends BaseMonteCarlo with NearestNeighbour with Mesure with ScatterPlot {

  var cluster = new Array[Boolean](vectors.size);
  var spins = new Array[Int](vectors.size);
  val n = vectors.size;
  var rand = new Random();

  var vecs = vectors.map(x => (knnPlusDist(x, vectors, k)))
  var meanAvr: Double = 2 * math.pow((3 * minAverageDistance + maxAverageDistance) / 4, 2)

  //calculate all forces Jij and store it into the tupel
  vecs = vecs.map(x => x.map(y => { (calcJij(y._1 * y._1, meanAvr, k), y._2) }))

  /**
   * Start the methode
   * @param mcSteps compute mcSteps 
   * @param T 		at temperature T
   * @return
   */
  def runMethode(mcSteps: Int, T: Double): Double = {
    init()
    val thermSteps = mcSteps
    for (i <- 1 to thermSteps)
      monteCarloStep(T);
    for (i <- 1 to mcSteps) {
      monteCarloStep(T)
      mesure(spins, n, q);
    }
    //figure(2)
   // scatterPlot(vectors, spins)
    val moments = computeAverages()

    return (n / T) * (moments._2 - (moments._1 * moments._1))
  }
  /**
   * bring system in start 
   */
 override def init() {
    //initialize all spins (Hot start)
    super.init()
    spins = spins.map(x => rand.nextInt(q));
   // figure(1)
    //scatterPlot(vectors, spins)
  }

  /**
   * compute one montecarlo step
   * @param T temperature
   */
  def monteCarloStep(T: Double): Unit = {

    createCluster(T)

  }
  /**
   * grow the cluster from a random startpoint i
   * @param T
   */
  def createCluster(T: Double) = {
    cluster = cluster.map(x => false);
    var newSpin = rand.nextInt(q - 1);
    val i = rand.nextInt(n);
    if (spins(i) >= newSpin)
      newSpin += 1;
    growCluster(i, newSpin, T)
  }

  /**
   * recursive function to try add a point out of the nearest neighbors
   * @param i
   * @param newSpin
   * @param T
   */
  def growCluster(i: Int, newSpin: Int, T: Double): Unit = {
    cluster(i) = true
    vecs(i).foreach(x => tryAdd(i, x._2, x._1, newSpin, T))
    spins(i) = newSpin
  }

  /**
   * Try to add a new poit to the cluster
   * @param i
   * @param j
   * @param Jij
   * @param newSpin
   * @param T
   */
  def tryAdd(i: Int, j: Int, Jij: Double, newSpin: Int, T: Double): Unit = {
    if (cluster(j) == false && spins(i) == spins(j)) {
      val pij = 1 - math.exp(-Jij / T); //grow propability
      if (rand.nextDouble() < pij)
        growCluster(j, newSpin, T); //grow the cluster from point j
    }

  }
  /**
   * Methode to calculate the strenght between the nodes i and j 
   */
  def calcJij(dist: Double, avr: Double, k: Int): Double = 1.0 / k * math.exp(-dist / avr);

}